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ABSTRACT 

The  last  few  years  have  seen  a  rapid  development  of  the  so-called  Lattice 
Algorithms  for  the  fast  solution  of  finite  data  least  squares  problems.  While 
a  fast  algorithm  has  been  given  for  finite  data  Instrumental  Variable 
Recursions,  as  yet  no  finite  data  lattice  schemes  have  been  given.  In  this 
work  a  lattice  algorithm  is  derived  for  a  finite  data  Instrumental  Variable 
Recursion  and  its  use  in  both  ARMA  and  ARMAX  Time  Series  models  is  indicated. 
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SIGNIFICANCE  AND  EXPLANATION 


In  the  field  of  signal  processing  in  Electrical  Engineering  one  of  the 
basic  problems  is  the  design  of  computationally  fast  algorithms  for  the 
solution  of  least  squares  problems.  For  example  in  the  problem  of  adaptive 
equilization  for  a  communication  channel.  A  signal  is  transmitted  over  a 
communication  channel  and  received  in  distorted  form.  It  is  then  passed 
through  a  linear  filter  called  an  equalizer  that  estimates  the  original 
transmitted  signal.  Usually  the  one  equalizer  must  deal  with  several 
communications  channels  and  should  have  different  parameters  for  each  one. 
Since  the  parameter  alteration  must  be  done  in  real  time  the  need  for  fast 
algorithms  is  apparent.  The  equalizer  problem  can  be  posed  as  a  least  squares 
problem.  Actually  fast  algorithms  for  least  squares  problems  have  been 
recently  given.  However  for  some  signal  processing  problems  (and  equalization 
problems)  other  estimation  schemes  must  be  used.  This  article  provides  a  fast 
algorithm  for  one  such  scheme,  an  Instrumental  Variable  scheme. 
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FINITE  DATA  LATTICE  ALGORITHMS 
FOR  INSTRUMENTAL  VARIABLE  RECURSIONS 


* 

by  Victor  Solo 

0.  Introduction 

A  fundamental  position  in  Signal  Processing  and  Time  Series  modelling  is 
occupied  by  the  Autoregressive  (AR)  model  or  all-pole  linear  model  driven  by 
white  noise.  The  infinite  data  least  squares  problem  (i.e.  estimation  of  the 
AR  parameters  by  minimizing  an  expected  squared  error)  yields  a  set  of 
Toeplitz  efjuations  for  which  the  Levinson-Whittle-Robinson  algorithm  provides 
a  fast  solution.  This  algorithm  can  also  be  re-expressed  as  a  lattice 
algorithm  (IxA):  in  this  case  it  is  Burg's  algorithm  (see  Hahhoul  [3]).  When 
only  finite  data  is  available  the  theoretical  autocorrelations  can  be  replaced 
by  estimated  ones  but  then  the  least  squares  nature  of  the  solution  is  lost. 
Recently  Lee  et.  al.  [41,  [51,  have  shovm  how  a  fast  algorithm  (a  lA)  may  be 
constructed  for  a  solution  to  the  finite  data  least  squares  problem  (i.e. 
minimization  of  a  sum  of  squared  errors):  for  a  readable  discussion  see  also 
Sartorlus  and  Shensa  [8]  and  Shensa  [9). 

For  infinite  data  Instrumental  Variable  (IV)  estimation  of  the  AR 
parameters  in  a  scalar  autoregressive  Moving  Average  (ARMA)  model,  Carayannis 
et.  al  [2]  have  provided  a  fast  algorithm  recently  -  it  is  an  extension  of  the 
one  of  Levinson.  For  finite  data  IV  parameter  estimation  in  general,  Ljung 
et.  al  [6]  have  recently  given  a  fast  algorithm  but  it  is  not  a  LA.  In  this 
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work  a  finite  data  lA  is  given  for  a  general  IV  problem  and  its  use  with  ARMA 
and  ARMAX  Time  series  models  is  indicated. 

Section  1  contains  a  derivation  of  an  infinite  data  LA  for  IV  estimation 
of  the  AR  parameters  in  a  block  identifiable  ARMA  model:  this  discussion 
which  is  much  simpler  than  the  finite  data  case  reveals  the  basic  idea  used  in 
all  the  derivations.  In  Section  II  an  infinite  data  LA  is  derived  for  a 
multivariate  model.  In  Section  III  a  general  infinite  data  IVLA  is 
presented.  In  Section  IV  it  is  shown  how  this  algorithm  may  be  applied  to  an 
ARMAX  (Autoregressive  Moving  Average  Exogenous  (i.e.  input))  Time  series 
model. 
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I.  The  ARMA  model;  infinite  data  IVLA 


Here  we  consider  the  estimation  of  the  AR  parameters  of  an  ARMA  Time 
Series  model  from  the  autocorrelations.  Actually  the  scalar  version  of  the 
lattice  algorithm  derived  in  this  section  is  implicit  in  Carayannis  et.  al  [2] 
article  which  itself  was  implicit  in  the  work  of  Risannen  17)  (see  C  below). 
However  the  present  derivation  is  especially  simple  and  provides  a  stochastic 
interpretation  for  various  quantities  appearing  in  the  algorithm.  Further, 
the  arguments  used  here  guide  the  more  complex  derivations  in  the  finite  data 
case  of  Section  11  below. 

A.  Preliminaries 

It  is  well  known  that  the  AR  parameters  of  a  multivariate  ARMA  (p,q} 
y^  (dimension  r)  obey  the  equations 

“  -R-q-,,-q.p  da) 

where,  with  »  E(^  Rp  nonsymmetric  Toeplitz  matrix 

R_  **  I  R  _  R  _  .  • .  •  R  _  I 

_p  — q  — q_i  — q-p+1 


^q+1 

-^+p-i 


...  R 


R  .  *  (R  .R  .^...R  ) 

— q-1,-q-p  — q-1 — q-2  — q-p 


(A?  ...  a'*’)  . 

—1  — p 


These  equations  may  be  alternately  written 


where 


A  H  »  -R 

-p  -p  — q-1  »-q-p 


(A  ...  A.)  »  A  J  ,  H  »  J  R 

-p  —1  -p-  -p - p 


0  ....  0  1 


1  I  .  .  .  0 
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Now  Up  is  a  block  Hankel  matrix.  If  =  rank(j^)  =  rank(R^)  =  pr  then 

the  ARMA  model  is  called  block  identifiable  (Akaike  [1])  and  the  p  -  AR 
matrices  may  be  obtained  by  solving  the  pr  equations  (1).  For  the  present 
discussion  the  assumption  of  block  identif lability  will  not  be  restrictive. 
Thus  it  is  henceforth  assumed  that  n  =  pr. 

n 

B.  Derivation 

To  derive  a  lattice  algorithm  begin  by  introducing  the  forwards 
prediction  error 


e  (n)  =  y  +  Z”  ^  y  ,  (2a) 

— m  ^  1  —a  -^-1 

where  =»  ...  is  defined  by  an  orthogonality  condition  that  ensures 


(when  map)  equations  (1)  hold:  namely 

E(e  (n)(y*  ,  ...  y* 

— m  -*n— q— 1  -^-q- 


.)) 


(2b) 


i,  • 


A™R  a  -R  . 

- ffl  — q«1 ,  -q-m 

It  will  be  convenient  to  introduce  the  IV  sequence  ^  *  in-q  that  (2b) 
can  be  written 


E(e(n)(z'...z'  ))=0  .  (2c) 

“Ui  —^••1  in 

It  will  be  seen  in  what  follows  that  we  can,  to  a  large  extent,  forget  that 
^  a  in-q'  will  also  prove  convenient  to  denote  the  fact  that  e,„(n)  is  a 
linear  combination  of  ^  ...  in-m  with  the  coefficient  on  being  by 

writing 


A.'"’  ■  4  ®  -  >'^‘4-1  •••  4-.’  • 

Then  clearly  ^(n)  is  determined  uniquely  by  the  two  requirements  (2c),  (3). 

The  next  step  in  a  lattice  derivation  is  to  find  an  order  update  equation 
for  ^(n).  For  this,  introduce  the  backward  prediction  error 

-  Zn-m-1  ^  ^1  Cl-t  ^-1 

where  B™  =  (B™  ...  b"*)  is  determined  by 

**“i  '  tn 
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E<r  (n-1)(z* 

— m  -n-1 


z*  ))  =  0 
— n— m  — 


(4b) 


that  is  b”  obeys 


or 


=  (B*"  J)R  =  -r 

—  — m  — - B  —  -q+ni»-<j+1 


B^R  =  -R 

- m  —  -q+l /-q+m 


(4c) 


(In  the  scalar  case  cf.  expression  (f)  of  Carayannis  et.  al.  [2]). 


Note  again  the  uniquely  determined  by  two  conditions  (4a), 


(4d) 


r  (n)  '  y  6  s 
— m  n,n+m-1 


Now  recall  that,  by  definition 


E(e  (n)(2' 
-tn+1  -n-1 


z*  z*  ,))  =  0 
—n-m  — n-m-1  — 


(4d) 


(5) 


Next  observe  that  for  any  K,  e^(n)  +  K  r^(n-T)-y„  6 
choose  to  satisfy 

0 

-m  — n-m-1 

.r 


E(e  (n)  +  K  r  (n-1 ))z*  , 

-m  - m  —n-m-1 


i  •  e  • 

where 


K  «  K 


—m  —in 


G  =  E(e  (n-1  )z'  , ) 

— m  m  — n-m— 1 

»  E(r  (n-1  )z'  ) 

“TO  “TO  1 

then  we  must  have  (via  (2c),  (3))  that 

e  (n)  =  e  (n)  +  r  (n-1) 

— m+1  — m  -m+1  — m 

This  is  the  desired  order  update  for  ej||(n).  (Note  the  notation 

Z~^  =  (r^)"\) 


(6a) 

(6b) 

(6c) 

(7a) 


Note  from  (2a),  (6b)  that 


G*^  =  R  ,  +  A™J  R 

— tn - q-m-1  - q-1,-q-n 


(in  the  scalar  case  cf.  expression  (9)  of  Crayannis  et.  al.  [2]). 

On  the  other  hand  comparing  coefficients  of  lagged  ^'s  in  (7)  yields 


\f1 


(7b) 
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A®'*'’  =a“  +  K^.  ^  1<j<m 

-0  -]  -m+l  -3 


(7c) 


Next  an  order  updating  formula  for  £_(n)  is  given.  Again  beginning  with 


E(r  ,{n)(z*  ...  z*  , ))  =  0 
— m+1  — n  — n-m-1  — 

we  see  that  for  any  ^ 

r  (n-l )  +  K  e  (n)-y  63 

-m  ~  n,n+m 

so  if  we  choose  to  obey 

E(r  (n-1 )  +  K  e  (n))z'  =0 


X .  e . 

where 


K  =  K®  =  -G®  E 
— m 

G®  *=  E(r  (n-1  )z  ' ) 


I  =  E(e  (n)z' ) 
— m  — m  — n 


we  must  have«  as  required 


r  (n)  *  r  (n-1)  +  K  e  (n)  . 

— m  — m  —m+l  — m 


Again  observe  that  from  (9a),  (9b) 


(8) 


(9a) 

(9b) 

(9c) 

(lOa) 


—  -q+m+1  *  •2-  —  -q+i,-q+m 

(in  the  scalar  case  cf.  expression  (11)  of  Carayannis  et»  al.  [2]).  On  the 
other  hand  comparing  once  more  coefficients  of  lagged  j^’s  in  (10a)  gives 

(10b) 
(10c) 


m+1  e 
,  =  K®  , 
—m+1  — ro+1 


„m+1  ,  „e  ;:in 

B.  =  B.  +  K  .A. 
—j  —j  -m+1  — 3 


An  order  update  can  now  be  obtained  for  E®,  2*^.  From  (7a),  (9c)  calculate 

— m  — m 


2"  , 
—m+1 


E(e  {n)z’ ) 
-m+1  — n 


E(e  (n)z')  +  .  E(r  (n-l)z') 

— m  — n  —m+1  ~m  — n 


2®  +  ,  G® 

— m  -m+1  -m 


Similarly 


2®  -  (I  -  )E®  by  (9a) 

—m+1  —  —m+1  —m+1  — m 


(I  -  K®..  • 


^+1  -  ^+1  ^+1 


(11) 


(12) 
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Initial  conditions  are 


So'"’  ■  in  ■  io‘"’ 


0  ]e 

Note  from  equations  (11  )f  (12)  that  in  the  scalar  case  Z  =  Z  , 

mm 

C.  Connections  with  Triangular  Decompositions. 

Recalling  equation  (4c)  consider  the  following  indent! ty 


L  ,  R  ,  V  , 
— m+l  — m+l  — in+1 


^+m,q+1 


R  R 

- q+m,-q+1  - q 


where 


D  .  =  R  +  B  R* 

— m+1 - q  —  -q+m,q+1 


R  N  =  -R’  ,  . 

-m  —  — q+m,q+1 

Clearly  recursive  application  of  this  identity  Iseginning  with  i  =  ill » 

_D^  =  will  produce  a  triangular  decomposition  of  ^  for  each  m.  Just 
such  a  scheme  was  given  by  Risannen  [7]  and  equations  (7),  (10),  (12)  form  one 
half  of  that  scheme.  The  other  half  consists  of  a  triple  set  to  produce  N^. 
These  equations  are  obtained  by  replacing  ^  in  (7),  (10),  (12)  by 
(see  Risannen  [7]). 

To  show  the  equivalence  of  (7),  (10),  (12)  to  Risannen 's  equations  make 

the  following  notation  changes: 

A.®A  .fB,  =B,  .  0<i<  m-1 

—1  -m-i  — i  ^,+1 

A  Rl 

then  (7b),  (7c)  read  (talce  A  =0  s  <  0) 

Am+1  -m+1  ,m  ^  3® 

A.  =»A  -  .=*A..+K  .. 

— i  m+1-i  — m+1-i  -m+1  -tn+l-i 


(13a) 


^+k’^,b"’,  0<i<m 

-^-1  —m+1  — m-i 
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II,  ARMA  model;  finite  data  tVLA 


To  use  the  scheme  of  Section  I  the  ^  must  be  replaced  by  estimates 
Then  however  the  resultiny  estimates  no  lonyer  obey  infinite  data 
orthogonality  conditions  such  as  (2b),  (4b).  Now  an  algorithm  is  given  which 
does  obey  finite  data  versions  of  these  orthogonality  conditions.  It  will  be 
called  a  finite  data  LA. 

The  argument  is  simplified  by  using  a  setting  such  as  that  of  Shensa 
[9].  Here  though  a  vector  version  of  that  discussion  is  needed.  If  is 

an  r-dimensional  discrete  time  series  with  ^  ^  t  <  0,  introduce  the 

infinite  data  matrix 


1  =  (Yo'  Xi'  l2  * 

Also  introduce  the  shift  operator  (i.e.  infinite  matrix) 

'i'’  y.  =  Xo'  -^1'  ••• 

so  C  <i.  =  (y^,  X3'  •••  • 

Note  that  ?  =  I  but  C  /  I.  Define  now  a  "matrix  inner-product" 

between  two  data  matrices 


<x,  y>  =  E"  X.  y! 


(14) 


and  observe  the  fundamental  relation 


<  5"  x,£>  =  X.  y.'  =  E,  X.  y! 

—  ^  n  0  -1-1  1  -1-1 

=  <X,  C  (£>„  ,  . 

—  n- 1 


1  I 

^0  ^+1 


(15a) 


Notice  that  need  not  be  of  the  same  dimension  for  the  definition  (14) 

to  ma)ce  sense.  In  this  respect  note  for  example  that  if  ^  is  an  arbitrary 
vector  (of  dimension  that  of  )  then 

a'<x,ip  =  e”  a*  X.  y .  =  <a' X,  (/>  .  (15b) 

Also  if  w  is  an  infinite  data  vector  w  =  (aw,  ow^  aw^  ...  .  (15c) 

It  will  be  convenient  to  denote  C  as  so  ij_  =  ij_q.  Thus  we  will  use 
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-Q 

the  IV  matrix  ^  =  C  K  ~  IL  _q’  Finally  the  following  time  update  formula 


will  be  useful 


(15d) 


A.  Lattice  Equations;  Statement 

Proceding  by  analogy  with  the  infinite  data  scheme  introduce  the  forward 
prediction  error 

e^(n)  (16) 

the  kth  component  for  ^(n)  is  ^(1t»n)  =  ^  A”(n)^  Note  that  the 

^(n)  are  r  *  r  matrices.  These  coefficients  are  determined  by  the 
orthogonality  conditions 


<e  (n),  2  .  >  =0  1  <  i  <  m 

— ffl  - in 

equivalently,  putting  A"'(n)  •  (A"*(n)  ...  A™{n)) 


where 


A“'(n)R  (n) 
—  — m 


R  ,  (n) 

- q-1 ,-q-m 


...  z  )' 

— k-m 


(17a) 


(17b) 


R  ,  (n)  =  E-  Xi,(5j,  1  •••  ^ 

- q-1,-q-m  0  ^  -k-1  —k-m 

Thus  equations  (17)  are  simply  the  equations  satisfied  by  the  mth  order, 

finite  data,  IV  estimate  of  the  AR  parameters.  It  is  assumed  that  R^(n)  is 

of  rank  rm  for  all  n:  If  E(R„(n))  is  of  rank  rra  then  this  event  occurs 

—in 

with  probability  one.  When  q  »  0  so  that  ^  ^  equations  (17)  are  the 

finite  data  least-squares  equations. 

Let  us  observe  that  ^(n)  is  determined  uniquely  by  two  conditions; 
(17a)  and  (17c) 

e^(n)  -  e  -  span(^  j/_^)  .  (17c) 
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This  is  so  since  equation  (17a)  (or  (I7b))  specifies  the  coefficients  of 


iL  -1  it  ~m  uniquely. 

Now  introduce  a  bacliwards  prediction  error 

r  (n)  =  w  -  e"*  b"  .in)u  . 

where  £^(n)  is  determined  uniquely  by 

r  (t\)  -  u  6  S 
— m  — m  O.m—1 

<r  (n),  Z  .>  =  0  0<i<  m-1 

— m  - in 


(18a) 


{18b) 


Now  the  lattice  equations  can  be  stated  as  follows:  A  derivation  is 
given  in  B  below.  Note  that  the  equations  are  easily  modified  if  the  inner 


product  is  changed  to 


rn  n-i 


where  0  <  w  <  1.  This  is  a  suitable  form  for  time  varying  equations. 
Finally  it  is  clear  by  inspection  that  the  algorithm  requires  0(m) 
operations. 


ORDER  UPDATES 


e  (n)  =  e  (n)  +  K*^(n)C  V  (n-1  ) 
— m+1  — n  — m  -m 


r  ^,(n)  =  C  ^r  (n-1 )  +  IC®(n)e  (n) 
— TO+1  — m  — m  — m 


(19a) 


(19b) 


K’^(n)  =  -G’^(n  )!“*■(  n-1) 
— m  — m  — n 


K®(n)  =  -G®(n)E"®(n) 

-HR 


E®(n)  =  (I  -  K'^(n)l^(n))E®  ,(n) 


£*^(0)  =  (I  -  K®(n)^(n))E’^  ,(n) 

—TO  — “  —TO  —TO  —in*  1 


Initial  conditions  for  these  are 


^(n)  =  Ej(n)  -  =^(n-1)  *  ^ 


(19c) 


(19d) 


(19e) 


(19f ) 


■ 
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TIME  UPDATES: 


m  >  1 


first 

wo  calculations  must  be  performed  before  the  next  two: 

G'^{n) 

»  G’^(n-I)  +  e  {n,n)&^  ,(n-1)(1  -  C  ,(n-1,n-1))  ' 

“in  “in  -^n— 1  ro—i 

(20a) 

— m 

G®(n) 

=  G®(n-1)  +  r  (n-1,n-1)6®  ,(n)(1  -  C  ,(n-1,n-1))  ^ 

—m  “in  “11“  1  m— 1 

(20b) 

“in 

6®  (n) 

=  6®'  (n)  -  6*^'  (n-1)I~'^(n-l)G®(n) 

— m-i  -111-1  —III  -HI 

(20c) 

-m 

—TO 

=  'S*^',(n-1)  -  6®'  (n)r"®(n)G’^(n)  . 

— m-1  — tn-1  — Bi  — m 

{20d) 

The  initial  conditions  are: 


^(n)  »  ^(n-1 )  +  y  ,  z’ 

— 0  -*«-1  — n 

(20e) 

oJ(„) 

(20f ) 

(20g) 

-z;.,  . 

(20h) 

Initializing  values 

for  the  whole  recursion  are 

lj(0)  -^(0)  =^z^ 

^(0)  =  0  =  G^(0) 

—0  —  *“0 

«J(0)  =  0  ,  6j{0)  . 

Finally  updates  for 

the  scalar  C  (n,n)  are 
m 

^m+1 

(n,n)  =  £  (n,n)  +  (n,n) 

m  —m  — m+i  — m+i 

(21a) 

CpCnjn)  »  •  (21b) 
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B.  Derivation 


(1)  Order  updates  for  ^(n)# 

To  find  an  order  update  equation  for  begin  by  observing 


-1 


<?  r  (n-n,2  •>  =  <r  (n-1),i;  z  .  , 

— m  ““  “1  n  “1  I 


=  <r  (n-1),z  > 

— in  —  -1+1  n-1 


=0  1<i<m  by{18b) 


(22) 


In  view  of  (17a),  (22)  consider  that  for  any 

S.<"'  -lo® 

<e  (n)  +  K  5“’r  (n-1),  z.>=0  1<i<m. 

— m  —  —in  —  —in  — 

Thus  if  we  choose  JC  so  that 

-1 


<e  (n)  +  K  C  r  (n-1 ),  Z  „  *  0 

— tn  —  ~tn  —  -ni-1  n  — 


i  •  6  • 


K  =  -G'^(n)E’'^(n-1 ) 


(23a) 


where  G'^(n)  =  <e  (n),  z  ,>  (23b) 

-tn  — m  - m-1  n 

E'^(n-I)  =>  <t"’r  (n-1),  z  „  ,>  =  <r  (n-1),  z  >  (23c) 

“  “iR"i  n  n~i 

then  we  must  have  (19a).  Note  that  when  using  this  update  only  the  nth 
component  need  be  calculated.  Similarly  we  can  find  equation  (19b)  with 


K®(n) 

-in 

=  -G®(n)E"®(n) 

(24a) 

G®(n) 

“W 

“O  n 

(24b) 

E®(n) 

=  <e  (n),  Z->^  . 

—in  —0  n 

(24c) 
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(2)  Order  Updates  for  E®(n),  25(n) 

■  .  ■  M,  .  .—I  I-.-  ■  —  ■ 

First  apply  (19a)  to 


E'en)  -  <a  (n),  2-> 
—in  —n  -K)  n 


<e  ,{n),  ,(n)  <  C~’r  ,(n-1),  z^> 

-m-l  —0  n  -ra-1  -m-1  -O  n 


E®  ,(n)  +  ,(n)G®  An) 

-0-1  -0-1  — m-1 


»  (I  -  ,(n)K®  ,(n))E®  (n) 

—  — m-1  — m-1  -m-1 

which  is  equation  (19c).  Similarly  {19d)  follows. 

(3)  Time  Updates  for  ^(n),  r,^(n) 

To  obtain  time  updates  for  G®(n),  G*^(n)  (the  usual  next  step)  it  will 

— ttt 

be  necessary  to  find  time  updates  for  should  be  noted  that 

r  ®  * 

in  the  least  squares  case  (q  »  0),  Z.~  ^  that  ^  ’W>is  is 

easily  seen  by  onserving  that  ^  *  i£. -m-i  (23b)  can  be  replaced  by 

C"’r  (n-1)  since  c"’r  (n-l )  -  c"V  «  S,  Thus 

G’^(n)  »  <e  (n),  c”^r  (n-1)>  .  Similarly  ^(n)  *  <5  (n-1),  e_(n)>^. 

— m  — m  — m  n  ~m  — m  — m  *» 

To  find  time  updates  for  we  are  guided  by  conditions 

(17a),  (17c).  Consider  then,  using  (I5d),  that 

<e  (n-1),  z  ,>  »  <e  (n-1),  z  <  >„  “  «  (n,n-1  )^  .  .  (25) 

—  -1  n-1  — m  ——in  -^n  -nt— i 


This  leads  us  to  introduce  the  infinite  data  vector 


(26a) 


with  c'"’’’(n)  »  (c““’(n)  ...  c®“,  (n))  (and  hence  I  ,(n))  determined  by 
—  —0  —m-1  m-1 

<5  An),  Z  ,>  »  2'  0  <  i  <  m-1  .  (26b) 

m-1  - in  -n-i 

Note  there  are  mr  unltnowns  in  c^“^(n)  determined  by  mr  equations 
(26b).  Now  observe  that 

-2-  -i^n  “  ^  -i+lVi 


^-1-(l-l)  *  ^-i 


1  <  i  <  m-1  . 
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Then  consider  that  (recall  (15c)) 


e  (n-1 )  -  e  (n,n-l)C~'c  ,(n-1)  -  (/„  6  S, 

— m  -m  m-l  -H)  1,m 

while  from  (25)  and  the  definition  of  e„(n-1 ) 

— n 

<e  (n-1)  -  e  (n,n-1  )4  ,(n-1),z  .> 

— m  — m  m-1  --in 

=  <e  (n-1 ),z  . >  =0  ,  ^  ^ 

— m  - in-1  —  1<i<m  . 


Thus  we  must  have  from  (17a),  (I7c) 

e  (n)  =s  e  (n-1)  -  e  (n,n-1)C  ,(n-1)  .  (27a) 

— ffl  — m  — m  m-1 

Similarly 

r  (n)  =  r  (n-1)  -  r  (n,n-1 )C  ,(n)  .  (27b) 

To  proceed  further  it  is  necessary  to  provide  both  order  and  time  updates  for 


C^(n).  In  the  least  squares  case  (q  =  0)  only  time  updates  were  needed. 
(4),  Time  tind  Order  Updates  for 

As  usual  we  are  guided  by  the  conditions  (26).  The  order  update  is 

considered  first.  Begin  by  recalling  three  facts 

<K  ,(n),  Z  .>  =  z'  .  0  <  i  4  m+1 

m+1  - in  — n-1 


<K  (n),  z  .  >  =  z’  . 

m  - in  —n-1 


0  4  i  4  m 


<r  (n),z  .>  =  0  0<i4m  . 

-m+1  - in  — 

Further,  for  any  r- vector 

C  (n)  +  n'  r., (n)  6  S  . 
m  —  —m+1  0,m 

Thus  if  we  choose  to  satisfy  the  r  equations 

<£  (n)  +  n'  r  .,(n),  Z  ^  ,>  =  z'  , 
m  —  —m+1  - m-1  n  — n-m-1 


n'  =  (n)E  ^  (n) 
— tn+1 


i.e. 

where  fi*^(n)  =  z  ,  -  <z  .,5  (n)> 

— m  —n-m-1  - m-1  m  n 

then  we  must  have  the  desired  order  update 

5  (n)  -  5„(n)  {n)r  (n)  . 

in*ri  in  —in  “mTi  -m+i 


(28a) 


(29) 


The  nth  component  of  this  equation  is  (21a). 
To  find  (21b)  tcdce  m  =  1  in  (26b)  so 


<5-(n),  z>  =  s' 
0  —  11  -11 


or 


so 


0*  0* 

(n)  <u,  z>  =  (n) <u,  z>  =  z' 

—0  n— 0  —  n  -n 

0  *  — s 

C-  (n)  =  z*  £-  (n)  whereupon 
-K)  — n  — 0 


the  first  component  of  this  is  (21b). 

We  turn  now  to  the  time  update.  This  time  consider  the  three  facts 

0  <  i  <  m 


(30) 


<C  (n),  z  . >  “2*4 

m  - in  -fi-i 


(n-1 ),  > 

m-l  —  -in 


=  z 


— n-i 


1  <  i  <  m 


<e(n),z  ■>  *0  1<i<m  . 

“HI  ——in  ““ 

Further,  for  any  r-vector  n_ 

m-1  —  — TB  0,m 

thus  if  we  choose  to  satisfy  the  r~equations 

,(n-1)  +  (n)»  Z/>>„  = 

©••1  —  — in  "K)  n  -n 


i«e« 


where 


n*  =■  &  ,(n)E^(n) 

—  — 1»-1  —in 

(n)  =  z  - 

— m-1  — n  — o  m-1  n 


then  we  must  have  from  (26)  the  desired  time  update 

-1, 


C  (n)  *  i  C  , (n-1 )  +6  ,(n)S  (n)e  (n) 

tn  w— 1  *^—1  —TO  *^n 


(28b) 


(31) 


Actually  there  is  also  an  order  change  in  this  equation  but  that  will  not 

cause  a  problem.  Now  it  is  possible  to  derive  an  interconnected  set  of  time 

0  r  6  r 

and  order  updates  for  the  G  (n),  G  (n),  6  (n),  6  (n). 

— m  — m  -m  -m 

(5).  Time  and  Order  Updates  for  G®(n),  G’^(n),  6®(n),  6’^(n) 

■ 

First  time  updates  are  given  for  G'^(n),  G®(n)  and  then  time/order 

“in  “m 

IT  0  r 

updates  are  given  for  6  (n),  5  (n).  Begin  by  applying  (27a)  to  G  (n). 
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G^(n)  =  <e  (n),  z  > 
-m  -m  - m-1  n 


<e  (n-Ui  z  >  -  e  (n,n-1)  <  (  5  ,  (n-1 ),  z  ,> 

— m  - m-1  n  — m  m-1  - m-1  n 


=  <e  <n-l  ),  z  ,  >  , 

— m  - m-1  n-1 


+  e  (n,n-1)(z'  ,  -  <5  ,(n-1),  z  >  ,) 

— m  — n-m-1  m-1  - m  n-1 


=  G'^(n-I)  +  e  (n/n-IjS*^  ,(n-1)  by  (28a). 


To  reproduce  equation  (20a)  take  the  nth  component  of  (27a)  to  see 

e  (n,n)  =  e  (n,n-1)(1  -  i  (n-1, n-1)) 

—m  “m  m— I 

and  substitute  this  into  (32):  equation  (20a)  results. 


Next  apply  (27b)  to  G  (n) 

— m 

G®(n)  =  <C  ^£_(n-1),  z^> 
-in  — O  n 


by  (24b) 


<r'r^(n-2),  z^>^  -  r^(n-1,n-2)  <  (n-l  ) ,  2_^>^ 


<C  r  (n-2),  z^>  . 

-in  —0  n-1 


+  r  (n-1,n-2)(2;  -  <C  ,(n-1),  z  >  ) 

^n  — n  m-1  — O  n 

+  r  (n-1  ,n-2)5*  (n)  by  (28b)  .  (33) 

—in  — m  — w-i 

Again  (20d)  is  found  from  the  nth  component  of  (27b)  namely 

r  (n,n)  ■  r  (n,n-1)(1  -  C  ,(n,n)) 

— m  — m  m-1 

and  (putting  n  to  n-1)  substituting  this  in  (33), 

Now  equations  (20a),  (20b)  must  be  completed  by  time /order  updates  for 

^(n),  ^(n).  Consider  then  applying  (31)  to  (28a) 

—in  — m  .  .r  .  ^ 

(n)  =  2*  -  <C  (n),  z  > 

— ffl  — n-m-1  m  - ra-1  n 


=  2'  1  -  5  ,  (n-1 ),  2  > 

~n-m-1  m-1  - m-1  n 


-m-l  — m  — m  - m-1  n 


X 


— n— B-i  m-i  —  “IB  n“i  — m— '  — n  — m 


6]['  (n-1)  -  (n)E"®{n)G'^(n) 

—m“i  -^R  — fli 


which  is  equation  (20d).  For  ^®(n)  using  (29)  in  (28b)  gives 


«®(n)  =  2*  -  <C”^C  (n-1),  2^> 
in  O  n 


2’  -  <?  K  ,  (n-1 ), 

— n  m-1  —0  n 


“  .  (n“1)E^(n-1  )<c“’r  (n-1 ),  z  > 

“fli—  1  —in  — o  n 


«®',(n)  -  6'^V(n-l)E"*‘(n-1  )G®(n) 
— m-l  — ffl-l  —in  -m 


which  is  equation  (20c). 

Finally  we  can  calculate  initial  conditions.  Now 

-  <V"’' 

^'(n)  »  ^  -  <r'yn-1),  Z^>^ 

however,  from  equation  (30) 

i-i>n 


2' .  - 

— n-1  -n  — O  —0 


which  is  (20g).  Similarly  (20h)  follows. 
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III.  General  IV  finite  data  LA 


A.  Statement 


Let  Yr.*  'o®  r-dimensional  discrete  Time  Series  and  consider  the 

— n  n  — n 

estimation  of  the  parameters  ^  =  (A^  ...  in  the  lagged  regression  model 

(or  finite  impulse  response  model) 

X  =J^fA.  y  .+n  n>0 

-n  1  —1  ^-1  — n 

by  means  of  an  IV  sequence  which  is  correlated  with  the  y.  but 

— n 

uncorrelated  with  the  coloured  noise  The  finite  data  IV  estimator 

A^(n)  is  obtained  by  solving  the  normal  equations 


)  =.  AP(n)E" 
)  —X  0 


4-1 


-k-p 


.4-r 


Note  i<0.  Now  a  lattice  scheme  for  solving  these 

equations  will  be  developed. 

Introduce,  as  expected,  the  forward  prediction  error 

(34) 


V  (n)  -  E®  A®  (n)u  . 

—TO  i  — -1 


Now  the  mr  parameters  ^(n)  *  (A®. (n)  ...  A®  (n))  are  uniquely  determined 

■^c  —^1  —xm 

2 

by  the  mr  orthogonality  conditions 

1  <  i  <  m  .  (35a) 

Thus  the  ^(n)  are  uniquely  determined  by  the  conditions  (35a),  (35b) 

(35b) 

Now  in  searching  for  an  order  update  equation  for  v  (n)  it  becomes 


<v  (n),  z  . >  *  0 
— m  - in  — 


V  (n)  -x.es, 

“in  1  ,m 


immediately  clear  that  we  must  also  use  the  bacltwards  prediction  error 

-  a  Bj(n)i, 

which  is  determined  as  in  Section  II  together  with 

e  (n)  ■  ^  ^{n)y  . 

“in  -K)  1  —  -t 
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t 


Now  the  finite  data  lattice  IV  recursion  can  be  stated 


Order  update;  Equations  (19),  ( 36 ) 


-1 


V  An)  =  V  (n)  +  K  (n)5  r  (n-1 ) 
— m+l  — m  — xm  -m 


K  (n)  =  -G  (n)Il"’^(n-1 ) 
— xm  — xm  — m 


(36a) 


(36b) 


Time  update  (m  >  1):  Equations  (20),  (37) 


K  (n)  =  K  (n-1)  +  V  (n,n)6'^  ,(n-1)(1  -  £  ,(n-1,n- 

— xm  — xm  — m  -m-l  m-1 


— xm  — xm  — m 

The  initial  conditions  are 


1))"’  .  (37a) 


K  -(n)  =  K  .(n-1)  +  z* 

— xO  — xO  -«  —n-1 


(37b) 


B.  Derivation 


To  find  the  order  update  for  v  (n)  it  is  only  necessary  to  repeat  the 
type  of  argument  used  in  Section  II  to  find  equations  (36)  where 


G  (n)  =  <v  (n),  z  _ 

— xm  — m  - m-1  n 


(38) 


To  complete  the  algorithm  it  is  sufficient  to  add  a  time  update  for 

As  in  Section  II  this  requires  a  time  update  for  V  (n)  which  is  easily  seen 

to  be  (cf.  derivation  of  (27a)) 

(39) 

n  m- 1 

Now  we  find  much  as  before  (cf,  (32)) 

(40) 

ID  -nn- 1 

However  the  nth  component  of  ( 39)  gives 

V  (n,n)  »  V  (n,n-1)(1  -  5  (n-1, n-1)) 

— m  —TO  m— 1 

and  substituting  this  in  (40)  gives  {37a).  Equation  (37b)  follows  from  (38) 
and  ( 1 5d ) . 


V  (n)  =  V  (n-1)  +  V  (n,n-1)C  ,(n-1) 

“W  -HI  lft—1 


G  (n)  *  G  (n-1)  +  V  (n,n-l)6*^  ,(n-1) 
— xm  -xm  -m  -m-l 
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IV.  Application  to  IV  estimation  of  scalar  ARMAX  Time  Series  Models 

Consider  the  estimation  of  the  parameters  a_,  _b  in  the  ARMAX  Time  Series 

model 


(1  +  a(L))x  =  b(L)u  +  (1  +  c(I.))e„  . 

n  n  n 

This  model  is  basic  in  Identification,  Econometrics  and  Time  Series  (L  is  the 
backshift  operator;  a(L)  =  a^L^;  ^  =  (a^  ...  a^)  etc.).  If  we  write  this 

model  in  regression  form  as 


X 

n 


=  Zf  a*. 


1  — xi  -^-i 


+  n 

n 


where 


(x  u  )•;  a  .  =  (a.b.)';  n  =  (1  +  c(L))e 
n  n  XI  r  1  n  n 


then  the  simplest  IV  scheme  uses  the  IV  sequence 

~  ^-p-1 

To  find  a  finite  data  lattice  recursion  we  begin  as  earlier  by 

Introducing  the  forward  prediction  error 

e  (n)  o™.  (n)j/  . 

xm  0  1  -XI  -i 

where  K,y  are  the  infinite  data  quantities 

X  =  Uq,x^,  ... 

y  =  * 

The  parameters  oi"*(n)  »  (a'**  (n)  ...  a”*  (n))  are  determined  by  the  equations 

I 


(41  ) 


<e  (n),  Z  ,  >  =0 

xm  - in 

where  Z  is  the  infinite  data  matrix 


1  4  i  <  m 


(42) 


z_  =  (£^,  •  • 

As  before  an  order  update  for  e„_(n)  necessitates  the  introduction  of 

xw 

the  backward  prediction  order 

£,(„)  .  (43) 

Now  it  is  clear  from  Section  II  that  calculation  of  this  requires  parallel 
calculation  (via  equations  (19),  (20))  of  the  forward  prediction  error 
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*  iio  -  ^ ‘"’ii  -i  • 

Now  recall  that  is  defined  by  the  equations 

<e  (n),  Z.>*0  1<i<m.  (17a) 

-m  - in  — 

However  since  the  first  row  of  y_Q  is  Xq  the  first  row  of  equation  (I7a)  is 
just  equation  (42).  To  put  it  another  way,  the  first  row  of  ^(n)  is 
Cj^jijCn).  Thus  the  solution  of  the  present  problem  is  necessarily  obtained  by 
the  algorithm  of  Section  II.  This  means  in  effect  that  to  estimate  ^  by  a 
lattice  scheme  we  must  fit  the  model 


y  =  A?  y  .  +  h 

i.e.  the  scalar  ARMAX  problem  is  naturally  imbedded  in  an  ARMA  problem  (of 

dimension  2).  In  the  least  squares  case  (^  *  this  idea  has  been  used 

by  Lee  et.  al.  [4]  though  they  did  not  point  out  its  natural  appearance.  Thus 

it  is  not  necessary  as  suggested  in  [4]  that  u^  obey  an  ARMA  model  for  the 

lattice  algorithm  to  work  (in  the  IV  case  or  AR  in  the  least  squares  case). 

All  that  is  needed  is  that  the  algorithm  produce  a  solution  to  equation  (42); 

which  it  does.  The  fact  that  it  does  so  by  fitting  an  ARMA  model  to  u^  (or 

AR  model  in  the  least  squares  case)  is  (in  this  respect)  incidental. 

Finally  it  is  worth  observing  that  the  discussion  of  this  section  will 

carry  over  to  the  case  of  a  multivariate  ARMAX  model 

(I  +  A(L))x  *  B(L)U  +  (1  +  c(L))6 
—  — n  —  — n  —  — n 

(here  A(L)  ®  A.  etc.)  in  which  the  block  Hankel  matrix  based  on  the 

—  I  —1 

2r-vector  sequence  ^  »  (j^  rank  2pr. 


Conclusion 

This  article  has  provided  a  finite  data  lattice  algorithm  for  the  fast 
solution  of  equations  arising  for  instrumental  variable  estimation  of 
parameters  in  Time  Series  models.  Two  basic  finite  data  algorithms  have  been 
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presented.  One  is  for  the  fast  solution  of  the  equations  satisfied  by  the  AR 
parameters  in  a  block  identifiable  multivariate  ARMA  model.  The  other 
algorithm  is  for  the  instrumental  variable  estimation  of  parameters  in  a 
multivariate  logged  regression  model  (or  finite  impulse  response  model)  with 
coloured  disturbance  noise.  It  turns  out  that  the  instrumental  variable 
estimation  of  the  transfer  function  parameters  in  scalar  ARMAX  Time  Series 
model  can  be  naturally  embedded  into  the  lagged  regression  model  so  the  second 
finite  data  lattice  algorithm  can  be  applied. 

A  natural  problem  for  future  investigation  is  the  derivation  of  a  finite 
data  lattice  algorithm  for  solving  equations  from  multivariate  ARMA  models 
that  are  not  block  identifiable.  This  would  include  the  general  ARMAX  case 
too. 


-23- 


REFERENCES 


[1]  H.  Akaike,  Markovian  Representation  of  Stochastic  Processes  and  its 

Application  to  the  analysis  of  Autoregressive  Moving  Average  Processes, 
Ann.  Inst.  ^tat.  Math.,  20,  p.363,  1974. 

[21  G.  Carayannis,  N.  Kalouptsidis,  D.  Manolakis,  Efficient  Algorithm  for 
Determination  of  A.R.  part  of  A.R.M.A.  model.  Elect.  Lett.,  656-657, 

1980. 

[3]  J.  Makhoul,  Stable  and  Efficient  Lattice  Methods  for  Linear  Prediction, 
l.E.E.E.  Trans.  Acoust.,  Speech,  Signal  Processing,  vol.  ASSP-25,  423- 
428,  1977. 

[4]  D.  T.  Lee,  B.  Friedlander  and  M.  Morf,  Recursive  Ladder  Algorithms  for 
ARHA  modelling,  Proc.  19th  IEEE  Conf.  Decision  and  Control,  1225-1231, 
Dec.  1980. 

[5]  D.  T.  Lee,  M.  Morf,  B.  Friedlander,  Recursive  Least  Squares  Ladder 
Estimation  Algorithms,  IEEE  Trans.  Acoust.,  Speech,  Signal  Processing, 
vol.  ASSP-29,  627-641,  1981. 

[6]  L.  Ljing,  H.  Morf,  D.  Falconer,  Fast  Calculation  of  Gain  Matrices  ^or 
Recursive  Estimation  Schemes,  Int.  J.  Control,  vol.  27,  1-19,  1978. 

[7]  J.  Risannen,  Algorithms  for  Triangular  Decomposition  c!  Block  Hankel  and 
Toeplitz  Matrices  with  Application  to  Factoring  Positive  Matrix 
Polynomials,  Math.  Comp.,  vol.  17,  147-154,  1973. 

[8]  E.  H.  Sartorius,  M.  J.  Shensa,  Recursive  Lattice  Filters  -  A  Brief 
Overview,  Proc.  19th  IEEE  Conf.  Decision  and  Control,  955-959,  Lee.  1980. 

[9]  H.  J.  Shensa,  Recursive  Least  Squares  Lattice  Algorithms  -  A  Geometrical 
Approach,  IEEE  Trans.  Autom.  Control,  vol.  AC-26,  695-702,  1981. 


VS/jvs 


security  Cl  ASSIFICAI  ion  of  This  PAGC  (■HAw  n«t»  Kntrml) 


REPORT  DOCUMENTATION  PAGE 

1.  hlpokt  numueh 

2.  GOVT  ACCESSION  NO. 

#2347 

4.  TITLE  (tuiJ  Subtitie) 

Finite  Data  Lattice  Algorithms  for 
Variable  Recursions 

Instrumental 

T.  AUTHORf.; 

Victor  Solo 

9.  PERFORMING  organization  NAME  AND  ADDRESS 

I  Mathematics  Research  Center,  University  of 

610  Walnut  Street 

Madison,  Wisconsin  '^3706 

Wisconsin 

11.  CONTROLLING  OFFICE  NAME  AND  ADDRESS 

U.  S.  Army  Research  Office 

P.O.  Box  12211 

■ 

Research  Triangle  Park,  North  Carolina  27709 

14.  MONITORING  VGENCY  NAME  ft  AOORESSCIF  dldcrant  From  Controtllng  Otiica) 

READ  INSTRUC  TIONS 
BEFOKK  COMPEETlNn  FORM 


3  RECIPIENT'S  CAT  ALOG  NUMQER 


S.  TYPE  OF  REPORT  ft  PERIOO  COVERED 

Summary  Report  -  no  specific 
reporting  period 


6.  performing  org.  report  NUMOER 


B.  contract  or  grant  NUM0ER('»; 

DAAG29-80-C-0041 


10.  PROGRAM  element,  project,  TASK 
AREA  B  WORK  UNIT  NUMBERS 

Work  Unit  Number  4  - 
Statistics  &  Probability 


>2.  REPORT  DATE 

March  1982 


u.  number  OF  pages 
24 


IS.  security  class.  (Oi  thit  report) 


UNCLASSIFIED 


is«.  declassification  downgrading 
schedule 


16-  OISTR'BUTION  statement  for  this  Report) 

Approved  for  public  release;  distribution  unlimited. 


t7.  OISTRIGUTION  statement  (ot  the  abetrect  entared  fr»  Bfocit  20,  i(  ditterant  /root  Report) 


19.  key  words  fConfinue  on  reveree  eida  H  necaaaary  and  idantUy  by  block  number) 

Time  series.  Signal  processing,  Instrumental  Variables,  Hankel  matrix, 
Becursive  least  squares,  ARMA  model.  Past  Algorithm,  Toeplitz  matrix 


20.  abstract  (Conitoua  oo  raver ea  eida  it  nacaaeaey  and  idertttty  by  bioek  rtianbar) 

The  last  few  years  have  seen  a  rapid  development  of  the  so-called  Lattice 
Algorithms  for  the  fast  solution  of  finite  data  least  squares  problems.  While 
a  fast  algorithm  has  been  given  for  finite  data  Instrumental  Variable  Recursions, 
as  yet  no  finite  data  lattice  schemes  hav^  been  given.  In  this  work  a  lattice 
algorithm  is  derived  for  a  finite  data  Instrumental  Variable  Recursion  and  its 
use  in  both  ARMA  and  ARMAX  Time  Series  models  is  indicated. 


OD  I  j An*7J  1473  edition  OF  I  NOV  ftft  IS  OBSOLETE 


UNCLASSIFIED _ 

security  classification  of  this  pace  Dmi*  Enl*r*d) 


